Clay hydroxyl isotopes show an enhanced hydrologic cycle during the Paleocene-Eocene Thermal Maximum

The Paleocene-Eocene Thermal Maximum (PETM) was an abrupt global warming event associated with a large injection of carbon into the ocean-atmosphere system, as evidenced by a diagnostic carbon isotope excursion (CIE). Evidence also suggests substantial hydrologic perturbations, but details have been hampered by a lack of appropriate proxies. To address this shortcoming, here we isolate and measure the isotopic composition of hydroxyl groups (OH−) in clay minerals from a highly expanded PETM section in the North Sea Basin, together with their bulk oxygen isotope composition. At this location, we show that hydroxyl O- and H-isotopes are less influenced than bulk values by clay compositional changes due to mixing and/or inherited signals and thus better track hydrologic variability. We find that clay OH− hydrogen-isotope values (δ2HOH) decrease slowly prior to the PETM and then abruptly by ∼8‰ at the CIE onset. Coincident with an increase in relative kaolinite content, this indicates increased rainfall and weathering and implies an enhanced hydrologic cycle response to global warming, particularly during the early stages of the PETM. Subsequently, δ2HOH returns to pre-PETM values well before the end of the CIE, suggesting hydrologic changes in the North Sea were short-lived relative to carbon-cycle perturbations.

warming is predicted to exceed 2°C by the end of this century even under the stabilisation (RCP6.0) greenhouse gas scenario 15 . Such warming will likely increase hydrologic-cycle intensity, thus exacerbating environmental stress 16,17 . Current models suggest an enhanced hydrologic cycle with elevated precipitation driven by higher evaporation rates 18,19 , a pattern which may already be emerging. However, while spatial changes in total precipitation amount can be modelled accurately, climate models are not yet precise enough to predict spatial variations in the frequency, seasonality, intensity or type of precipitation 17,20 . Testing whether current climate models can simulate the observed hydrologic responses to temperature during the PETM will aid in assessing future hydrological projections in response to global warming.
Large variations in clay mineralogy are recorded in sediments deposited over the PETM in response to changing climate conditions. A major unresolved observation is a conspicuous increase in kaolinite deposition at many mid-to high-latitude sites (e.g. refs. [21][22][23][24][25][26]. Kaolinite is a clay mineral that is generally assumed to form predominantly by intense weathering in humid tropical climates [27][28][29] . The origin and significance of PETM kaolinite deposition is controversial because it could result from either: (i) increased chemical weathering and kaolinite formation under warmer climate and elevated year-round precipitation 23,[30][31][32] , or (ii) enhanced mobilisation, physical transport, and exhumation of previously weathered material (e.g. laterite soils) during strong seasonal rainfall events 26,28,33 . Whereas the former interpretation of increased chemical weathering is consistent with a direct negative silicate-weathering feedback to elevated temperatures 34 , the latter of increased physical weathering provides a less direct feedback and may indicate a decoupling between climate and weathering intensity. The two processes are not mutually exclusive and evidence exists that both may have occurred at the time of the PETM in different settings.
To provide additional constraints on the hydrologic changes across the PETM, we measure the isotopic composition of clay mineral hydroxyl groups from a PETM section in the Sele Formation of the North Sea (well 22/10a-4 (57°44′8.47″N; 1°50′26.59″E; up to ∼500 m paleo-water depth); Fig. 1; refs. 10,33,35). During the early Paleogene, the North Sea was a restricted marine basin bounded by Scotland to the north, Greenland to the west and the Fenno-Scandian Shield to the east. The basin contains a highly expanded Paleocene-Eocene transition sequence, uninterrupted but for minor erosion at the base of thin turbidite sandstones 10 . Terrigenous input was high during the Paleogene and thought to be derived primarily from the Scotland Faeroe-Shetland landmass 36 , making this an ideal location to study PETM variability in the hydrologic cycle at a mid-paleo-latitude location (57°N). Furthermore, by isolating and measuring the hydrogen and oxygen isotope composition of clay hydroxyl groups 37 (δ 2 H OH and δ 18 O OH ; Methods) in addition to bulk clay δ 18 O (δ 18 O bulk ; Methods), our approach aims to avoid potential biases by clay end-member mixing and/or inherited isotopic compositions 26 . We thus expect this record to robustly reflect PETM variability in the hydrologic cycle.

Results and discussion
We compare all clay isotope results relative to the record of bulk organic matter δ 13 C from the same section 10 and previously reported clay mineralogical data 33 (Methods; Fig. 2; Fig. 3). The CIE onset in this section is located between 2614.3 and 2613.5 m 10,33 , whereas the post-PETM δ 13 C recovery is not captured in the available δ 13 C data and likely occurs higher in the section.

Trends in hydroxyl hydrogen and oxygen isotopic composition
Measured δ 2 H OH values decrease gradually by 13.7 ± 0.6‰ over the 15 m leading up to the CIE onset, reaching −67.8 ± 0.3‰ VSMOW (Vienna Standard Mean Ocean Water) at 2613.5 m (Fig. 2b). The CIE onset at 2613.5 m is followed by an abrupt δ 2 H OH decrease of 8.1 ± 0.4‰; δ 2 H OH decreases further throughout the early CIE, culminating in a minimum value of −87.0 ± 0.2‰ VSMOW at 2611.6 m. This δ 2 H OH excursion occurs only between 2613.5 and 2608.7 m, with δ 2 H OH       In contrast to δ 2 H OH , δ 18 O OH changes appear to be more closely related to changes in clay composition and δ 13 C (Fig. 2e), but overall correlation with clay oxygen isotope values is weak (e.g. for kaolinite, the strongest correlation, r 2 = 0.46; p-value < 0.01; n = 22). At the base of the section, δ 18 O OH averages 3.2 ± 0.3‰ VSMOW but increases to 10.7 ± 0.2‰ VSMOW (omitting a one-point minimum of 1.38 ± 0.13‰ VSMOW at 2617.4 m, which corresponds to a kaolinite content maximum and illite-smectite (I-S) content minimum; Fig. 3a

Estimation of clay hydroxyl source water
Finally, we estimate the hydrogen-isotope composition of clay mineral (neo)formation source water (δ 2 H OH-SW ) using mineral-specific 2 H fractionation factors 39,40 combined with measured δ 2 H OH and clay mineral content trends (Methods; Fig. 3e). Resulting δ 2 H OH-SW closely tracks the measured δ 2 H OH value (r 2 = 0.90; p-value < 0.01; n = 22). Absolute values are less negative (ranging from −2.5 to −50.9‰ VSMOW) and exhibit greater variability (48.5‰) relative to δ 2 H OH (34.1 ± 0.3‰), as the most enriched δ 2 H OH values correspond to depths with higher I-S concentrations, while the most depleted δ 2 H OH values correspond to depths with higher kaolinite concentrations, despite kaolinite having a larger fractionation factor.

Implications of bulk and hydroxyl isotopes on clay origin
The origin of the pervasive increase in PETM kaolinite deposition is debated. Early studies interpreted this trend as reflecting contemporaneous PETM kaolinite formation, indicating increased pedogenesis and weathering under a warm, humid climate with an intensified hydrologic cycle 23,[30][31][32] . More recent studies considered the time required for regolith kaolinitization-with estimates on the order of one million years 28 -and suggested instead that this increase in kaolinite represents increased erosion and exhumation of previously deposited kaolinite that formed well before the PETM 26,33,41 . However, such long timescales refer to the formation of kaolinite from unweathered bedrock but do not necessarily apply to the kaolinization of pre-existing clays in soil profiles. The soil transformation of 2:1-type phyllosilicates (e.g. smectite and illite) into 1:1 types (e.g. kaolinite) has been found to occur on the order of 100 days under laboratory conditions at 150°C 42 , and significant kaolinization of smectites can occur over a 5-10 ka period under tropical weathering conditions 43 .
Clay oxygen isotopes provide information about possible mineral sources and their genesis. As such, bulk oxygen isotope measurements of PETM-age clays have been analysed previously; for example, from the Bass River section on the New Jersey margin, which similarly contains a mineral assemblage of smectite, illite, kaolinite and quartz 26 . In that study, a strong correlation was observed between smectite and δ 18 O bulk (r 2 = 0.96) and between kaolinite and δ 18 O bulk (r 2 = 0.88). This suggested that δ 18 O bulk variability resulted from compositional changes in the proportion of kaolinite and smectite-each with constant δ 18 O bulk values-rather than changes to the δ 18 O bulk of the individual clay minerals across the PETM. Applying a similar analysis to our δ 18 O bulk measurements, we find that only kaolinite shows a similarly strong correlation (r 2 = 0.83) with δ 18 O bulk (Table 1). This implies that kaolinite across the section has a fairly constant δ 18 O bulk value, suggesting that kaolinite δ 18 O bulk does not necessarily capture isotopic changes across the PETM; this could be the result of kaolinite erosion from a pre-existing source or from kaolinization of existing soil profiles (where δ 18 O bulk will only be partially affected by the addition of hydroxyl groups).
By contrast, the relative amounts of oxygen contributed by smectite, illite, quartz and chlorite components do not strongly correlate with δ 18 O bulk , δ 18 O OH or δ 2 H OH , with r 2 values never exceeding 0.52 (Table 1). Kaolinite content-which contributes up to 75% of the oxygen and hydrogen for the hydroxyl measurements-also shows poor correlation to δ 18 O OH and δ 2 H OH , suggesting that the variation in the hydroxyl isotopes over the section is not driven by changes in the clay-sized mineral proportions.
We thus do not observe strong clay compositional effects in the section studied here. In fact, the δ 2 H OH-SW variability increases when measured isotopic values are corrected for clay composition (i.e. by mineral-specific fractionation factors), suggesting a greater change in source water isotopic composition than is actually measured. Furthermore, the measured δ 2 H OH-SW decrease (48.5‰) is consistent with n-alkane δ 2 H data from the same region (the Normandy Vasterival section), where a decrease of 60‰ was measured across the early PETM 14 . The observed changes in δ 2 H OH are thus best explained by the formation and alteration of clays during the PETM, rather than by the erosion of pre-existing deposits. In order for kaolinite formation to capture changing source water isotopic composition across the PETM, it must have formed via a faster mechanism than regolith kaolinization (as regolith kaolinization is estimated to take longer than the duration of the PETM 28 ). We therefore suggest that kaolinite was formed by the transformation of 2:1type phyllosilicate clays (e.g. smectite and illite) into 1:1-type clays (e.g. kaolinite) within soil profiles. This kaolinization process of smectitic precursors leads to the addition of new hydroxyl groups-which will isotopically reflect meteoric water-and has been shown to occur under warm, humid conditions at rates sufficiently fast to capture changes in rainfall isotope composition across the PETM 44,45 .
If, as argued above, clay formation and alteration during the PETM is the major contributor to δ 2 H OH in this section, then the effect of temperature on δ 2 H OH must also be considered, as temperatures in the North Sea region increased over the PETM. Current models for clay hydroxyl group H-isotope fractionation factors ( 2H α OH ) predict values less than one 39,40 , meaning that clay hydroxyl group H-isotopes are depleted relative to the hydroxyl group source water. Less isotopic fractionation occurs at higher temperatures, so we would therefore expect temperature increases to result in higher δ 2 H OH values over the CIE in contrast to the observed lower values.
The δ 2 H OH signal may be somewhat biased towards marine values by post-depositional isotopic exchange; however, in contrast to exposed terrestrial clays, there is little evidence for hydrogen isotopic exchange in deep-sea buried clays over ∼3 Myr outside of the <0.1 μm fraction 46 , so significant hydroxyl exchange is unlikely. Contamination by organic matter may also affect the accuracy of δ 2 H OH measurements; however, we consider this effect to be relatively minor to the overall δ 2 H OH trend (see the Supplementary Information for a more complete discussion).

Paleohydrologic implications of δ 2 H OH
Taking into account the weak or insignificant correlations between measured isotope values and clay mineral proportions, and the negative relationship between temperature and 2H α OH , we interpret the hydroxyl isotope record as capturing changes in the isotopic composition of the source water. Although the decrease in δ 2 H OH at the PETM onset is the largest observed change, both the measured δ 2 H OH and the calculated δ 2 H OH-SW show a clear trend towards more negative values in the period leading up to the PETM onset. This trend towards more negative δ 2 H OH values precedes the CIE, suggesting that hydrologic changes during the PETM were not solely a result of changes in the carbon cycle. The pre-PETM decrease in δ 2 H OH may instead reflect increased terrestrial input to the North Sea, as sea level in the basin was lowered by ∼100 m during this time 47 . However, there is additional palynological evidence from a close-by North Sea well, where significant changes to vegetation were documented prior to the CIE, implying the pre-PETM decrease in δ 2 H OH does reflect hydrologic changes which predate the CIE 48 .
The 8‰ VSMOW decrease in δ 2 H OH at the onset of the PETM is interpreted to reflect a decrease in the δ 2 H OH of rainfall. We propose that the mechanism causing this δ 2 H OH decrease is the amount effect 49 -an observed negative correlation between monthly mean rainwater isotope composition and the total amount of precipitation, typically seen in the tropics at low elevation due to increased precipitation intensity in those regions. Increased PETM precipitation amount is further evidenced in this section by the concurrency of the δ 2 H OH excursion and the increased abundances of the low-salinitytolerant dinoflagellate cysts (and Apectodinium), implying that δ 2 H OH accurately reflects the δ 2 H of rainfall during the PETM.
The measured δ 2 H OH decrease during the PETM may reflect high intensity precipitation events where rainfall is depleted in 18 O and 2 H; for example, those resulting from tropical cyclone activity 50 . TEX 86derived sea-surface temperatures for the North Sea suggest a temperature increase of at least ∼10°C across the CIE onset with temperatures exceeding 30°C 51 . Such high temperatures exceed the tropical convective threshold 52 , and may have contributed to the high rainfall intensity. Our δ 2 H OH results are consistent with previous interpretations that global warming during the PETM caused an intensified hydrologic cycle and northward migration of storm tracks 10,11,41 .
The δ 2 H OH values return to their pre-PETM baseline before the CIE termination, suggesting that increased rainfall intensity around the North Sea was short-lived; this too is consistent with prior δ 2 H measurements from the region, where short-lived increases in seasonality have been inferred from n-alkane δ 2 H in the Normandy Vasterival section 14 .
Similar hydrological behaviour during the CIE is suggested by mineralogical data from the Fur section in Denmark, where the total clay fraction, clay mineral proportions and the clay chemical index of alteration return to pre-PETM values before the end of the CIE 53 . Lithium isotopic composition measurements from the Fur section similarly show an excursion at the PETM onset before beginning to return to pre-PETM isotopic values during the CIE 54 . However, relative changes in proxies around the CIE onset differ in the Fur section compared to ours-in particular, only the total clay fraction seems to increase before the CIE onset, and kaolinite in this section is mostly present after the CIE onset. The Fur section is less expanded over the start of the PETM (i.e. the Stolleklint clay layer), suggesting that early hydrological changes may not have been fully captured.
Also of note is that the 8 ‰ VSMOW decrease in δ 2 H OH at the PETM onset appears to precede the CIE in this section. This may offer further evidence that hydrologic changes during the PETM were not a result of carbon excursion. However, the CIE onset recorded in this section is complex and potentially influenced by multiple confounding factors, including variable mixing of different carbon sources at time of deposition 55 (where carbon sources could be either marine or terrestrial, and may or may not reflect newly sequestered carbon), and possible variability caused by a pre-onset carbon isotope excursion, which has been recorded in other sections 56 (although is not present in this section's measured record).

Future Direction
This study uses the recently developed differential thermal isotope analysis (DTIA) method to report clay δ 2 H OH across the PETM 37 . Clay hydroxyl isotopic composition is not a well-developed area of research; future work will benefit from better constraints on hydrogen isotopic fractionation factors for clay hydroxyl groups ( 2H α OH ). There currently exist a limited number of 2H α OH values in the literature 39,40 , and reported values assume that fractionation is not influenced by variations in chemical composition (i.e. the interlayer metal ions present) or mineralogical composition (i.e. the ratio of illite to smectite within illite-smectite). New 2H α OH estimates for a suite of clay minerals under environmentally relevant conditions are needed to more accurately correct for changing clay composition and estimates of paleoprecipitation δ 2 H variability. Methodological improvements for the removal of organic matter without alteration of the clay hydroxyl isotope composition, or better constraints for the influence of organic matter on DTIA δ 2 H OH measurements, will also improve the certainty of this proxy. Similar clay isotope measurements should be made in other PETM sections to determine whether the δ 2 H OH values reported here indicate pervasive hydroclimate changes at the PETM onset or were driven by local effects in the North Sea Basin. We emphasise the importance in future studies of carefully evaluating the effects of changing clay composition on the δ 18 O bulk , δ 18 O OH and δ 2 H OH before an environmental interpretation is proffered.
In conclusion: we have applied a recently developed method for measuring δ 2 H and δ 18 O of the hydroxyl group in clay minerals to reconstruct past hydrological changes during the PETM. We observe a large decrease in δ 2 H OH coinciding with the CIE, which we interpret to reflect increased precipitation-possibly resulting from increased tropical cyclone activity-in response to global warming during the PETM. Our interpretation is supported by concomitant increases in the abundances of low-salinity dinoflagellates 10 , indicating increased runoff to the North Sea Basin. Hydrologic changes are short-lived relative to the CIE, as evidenced by δ 2 H OH values and low-salinity-tolerant dinoflagellate abundances returning to pre-PETM values before recovery of the carbon cycle. Taken together, these results imply increased kaolinization, silicate-weathering intensity, and hence humidity, in response to warming events in a greenhouse world. Comparing these and future PETM clay δ 2 H OH measurements to predictions from isotopically enabled models offers the promise of quantitative reconstructions of paleohydrologic conditions that can enhance our understanding of hydrology and pedology in a warming world. Moreover, clay hydroxyl group isotopes are a largely untapped proxy, which can greatly improve our understanding of clay provenance and paleoclimate in the past, provided the compositional changes in clay mineralogy are fully characterised.

Isotope notation
All isotope ratios are reported in traditional "delta" notation; for example, written here for 2 H/ 1 H as where 2 R is the 2 H/ 1 H ratio and "standard" refers to Vienna Standard Mean Ocean Water (VSMOW) for δ 2 H and δ 18 O of hydration water and clay bulk oxygen or to Vienna Pee Dee Belemnite (VPDB) for δ 13 C for calcite. All results are reported in units of "per mil" (‰).

Samples
We measured 22 samples in this study consisting of separated <4 μm size fractions of sediment samples across the PETM from the Sele Formation, taken from well 22/10a-4 (57°44′8.47″N; 1°50′26.59″E; up to ∼500 m water depth) located in the central North Sea. Bulk organic matter δ 13 C has been reported elsewhere for the same section 10 . In the late Paleocene prior to the PETM, δ 13 C remains relatively constant at −25 ± 1‰ VPDB (Fig. 2a). CIE onset in this section is located between 2614.3 and 2613.5 m 10,33 , with a large δ 13 C decrease to −30 ± 1‰ VPDB at 2612.0 m. Carbon isotope values remain low (−30.7 ± 0.3‰) until at least 2605.0 m; post-PETM δ 13 C recovery has not been captured with the available data and likely occurs higher than measured in this section.

Sample preparation and measurement-hydroxyl isotopes
Details of the initial preparation and isolation of the <4 µm samples are given in ref. 33. Before beginning isotope measurements, a test sample was heated using a heating ramp of 5°C per minute to 1030°C to identify the temperature range of dehydroxylation (Fig. 4).
Before making the hydroxyl isotope measurements, we first treated the <4 µm samples with 2 ml of 2 M sodium acetate/acetic acid buffer solution to remove any traces of carbonate. Following acid treatment, the samples were dried overnight at 40°C and then measured with the Differential Thermal Isotope Analysis (DTIA) system 37 , using a heating programme optimised for maximum separation of the different forms of water contained in the clay minerals (Fig. 5). All heating ramps used were at a 40°C / min rate to give sharp water peaks. Sample sizes varied between 16 and 31 mg and were adjusted to give a hydroxyl peak height of 11,000-14,000 ppm water in the optical cavity of the Picarro 2130 analyzer. The δ 18 O and δ 2 H were calculated by integrating the H 2 O, δ 18 O and δ 2 H traces for the hydroxyl peak and correcting them for background and calibration of the Picarro cavity ringdown laser spectrometer.
Three to four DTIA measurements were taken at each depth (at 2610.82 m, only two measurements could be taken, due to having limited sample). Errors (1 σ) were calculated from the standard deviation of the measurements at each depth.

Sample analysis-hydroxyl isotopes
To calculate the fractionation factor corrected δ 2 H values from the hydroxyl group δ 2 H, which should represent the approximate source water hydrogen isotopic composition, a number of assumptions are required. Not all clays have hydroxyl fractionation factors reported in the literature, so we must assume that the hydroxyl water was derived from a mixture of illite-smectite and kaolinite only. We then calculate the proportion of hydroxyl water that each clay mineral would contribute to the total hydroxyl water, by comparing the proportions of each clay in the assemblage to the relative weight percentage of hydroxyl in the clay chemical formulae, assuming a typical clay elemental composition. For illite-smectite, we assume that changes to the percentage composition of illite, K 0.65 Al 2.65 Si 3.35 O 10 (OH) 2 , and montmorillonite, Na 0.1 Ca 0.23 Al 2 Si 4 O 10 (OH) 2 (H 2 O) 11.66 (assuming the montmorillonite is 36.1 % water by weight), do not affect the fractionation factor, but do change the proportion of hydroxyl H-isotopes contributed. The source water values were then estimated by applying the fractionation factors weighted by the water amounts contributed by each clay: where: α 2H mineral is the water-mineral hydrogen-isotope fractionation factor, and δ 2 H SW is the estimated source water value for δ 2 H. Hydrogenisotope fractionation factors at 30°C are taken from refs. 39, 40. for kaolinite and illite-smectite respectively ( Table 2). We also estimate δ 2 H SW assuming compositions of 100% kaolinite and 100% illitesmectite, to account for the spread of possible values for δ 2 H SW in order to estimate the errors for the calculations (where the largest error is caused by not having fractionation factors for an average of ∼25% of the hydroxyl isotope contribution).
No reliable values for hydroxyl oxygen isotope fractionations in clays exist in the literature. As such, it was not possible to carry out this analysis for hydroxyl oxygen. Fractionation factors for bulk oxygen  Following this, samples were dried overnight at 65°C in a vacuum oven, loaded into a 21-well laser fluorination tray and dried for a further 72 h at 65°C in a vacuum oven. This tray was then inserted into a laser fluorination chamber, heated to ∼65°C under vacuum for ∼6 h, and allowed to cool to room temperature before being reacted with 40 Torr fluorine gas (F 2 ) for 15 h to remove any remaining hydration water ("pre-fluorination"). Pre-fluorination product O 2 was then removed under vacuum. To analyse each sample, the chamber was again charged with 40 Torr F 2 gas and same material was reacted by use of a 50 W CO 2 laser (Teledyne) following ref. 57; laser power was initially low to prevent sputtering and was gradually increased to ∼10% until no visible sample residue remained. Fluorinated compounds were then removed by liquid nitrogen, residual F 2 gas was passivated by KBr salt held at 200°C, and O 2 gas was purified by passage through a 3 m long gas chromatography column packed with 5 Å molecular sieve as described in ref. 58.
Analyte O 2 gas was measured on a ThermoFischer Scientific MAT 253 isotope ratio mass spectrometer (IRMS) operated in dual inlet mode; for each sample, gas was analysed in 4-acquisition blocks of 10 cycles each and averaged. Resulting δ 18 O values were corrected to the 2-point calibration SMOW/SLAP scale following ref. 59. by directly fluorinating and analysing VSMOW-2 and SLAP-2 standard material using a CoF 3 reactor as described in ref. 60. Uncertainty was taken as the difference between replicate sample aliquots.

Data availability
The clay isotope data generated in this study, as well as source data for all figures, are provided in the Supplementary Data.

Code availability
A Mathematica script which can be used in the analysis of raw Picarro L2130i cavity ringdown spectrometer data is included in the Supplementary Software.  isothermal steps were used to increase peak separation. The second peak varies in size across the samples, being between 12 and 51% of the size of the third peak, and thus is too small for accurate isotopic measurements. Including the second peak as part of the hydroxyl peak in isotope analysis leads to only small changes in reported values, and does not change the overall trends in δ 18 O OH and δ 2 H OH .